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Abstract 


Three explicit multigrid methods, Ni's method, 
Jameson's f Ini te- volume method, and a finite- 
difference method based on Brandt's work, are 
described and compared for two model problems. All 
three methods use an explicit multistage Runge- 
Kutta scheme on the fine grid, and this scheme is 
also described. Convergence histories for Inviscid 
flow over a bump In a channel for the fine-grid 
scheme alone show that convergence rate Is propor- 
tional to Courant number and that implicit residual 
smoothing can significantly accelerate the scheme. 
Ni's method was slightly slower than the 
Implicitly-smoothed scheme alone. Brandt's and 
Jameson's methods are shown to be equivalent In 
form but differ In their node versus cell-centered 
implementations. They are about 8.5 times faster 
than Ni's method in terms of CPU time. Results for 
an oblique shock/boundary layer interaction problem 
verify the accuracy of the finite-difference code. 
All methods slowed considerably on the stretched 
viscous grid but Brandt's method was still 2.1 
times faster than Ni's method. 

I ntroduction 

Three explicit multigrid methods now being 
used for solution of the Euler equations and occa- 
sionally for the Navier-Stokes equations are dis- 
cussed in this paper. The first method was 
developed by Ni In 19811 to accelerate the con- 
vergence of his own fine-grid Euler scheme. Sub- 
sequent work by Johnson and Chlma?, b generalized 
the method to other finite-difference fine-grid 
schemes and to viscous flows for practical turbo- 
machinery problems. A variation of Ni's method has 
also been used by Hall for Inviscid flows over air- 
foils. ^ A second method was developed by Jameson 
and Baker in 1 983"? to accelerate the convergence 
of the finite- volume multi-stage Runge-Kutta 
schemes developed by Jameson, Schmidt, and 
Turkel. 8 This method has been used very suc- 
cessfully for Inviscid two-dimensional flows over 
airfoils and three-dimensional flows over aircraft 
configurations ,9 and for two-dimensional viscous 

flows .10 
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Multigrid methods were first developed for 
elliptic problems and have been analyzed in detail 
by Brandt. 11 For the third method considered 
here Chima and Schaffer used Brandt's approach to 
develop a finite-difference multigrid Euler and 
Navier-Stokes code. Many of Jameson's Ideas were 
used on the fine grid so that a direct comparison 
could be made between the basic Runge-Kutta scheme, 
the multigrid scheme, and Ni's scheme. In par- 
allel, Turkel has revised Jameson's f ini te- volume 
multigrid Euler codes ( FL052MG and FL053MG) 
allowing comparisons with this method as well. 


The Intent of this paper is to develop and 
compare the Ni, Jameson, and Brandt types of mul- 
tigrid schemes. All of these schemes use the 
Runge-Kutta method as the basic algorithm. There- 
fore we will discuss the Runge-Kutta method, its 
associated boundary conditions and artificial vis- 
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volume forms. We will also discuss other con- 
vergence acceleration techniques Including spa- 
tially-varying time steps and residual smoothing. 
Details of the three multigrid methods will then be 
presented. Finally results for several inviscid 
and viscous flows will be used to demonstrate the 
relative effectiveness of the schemes and to point 
out areas where further work is needed. 


Governing Equations 

The two-dimensional unsteady thin-layer 
Navier-Stokes equations may be written in fully 
conservative form for an arbitrary coordinate sys- 
tem as follows: 


a t q = -J[3^E + a n (F - Re 1 S) - D] (1) 


where 
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e - p [C y T +- 2 (u f v )] is the total energy per 
unit volume. 
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p s (y - 1) ~ p(u 2 * V 2 )] is the static 

pressure and D is an artificial dissipation term 
to be described later. 


The viscous flux term S is given by: 
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In Eqs. (1) to (3) J is the transformation 
Jacobian: 


(3) 


For stability and shock capturing both schemes 
require an artificial viscosity term which will be 
discussed. Also both schemes use a local time step 
and implicit residual smoothing to accelerate con- 
vergence at the expense of losing time accuracy. 

In addition the finite-volume code uses enthalpy 

damping. 8 

Boundary Conditions - Finite-Difference Code 

For subsonic flow, inlet values of P 0 , T 0 , 

and v are specified. The upstream-running Riemann 
invariant R- = u-2c/(y-l) is extrapolated from 
the Interior and isentropic relations are used to 
compute the necessary flow quantities. For super- 
sonic flow all inflow quantities are specified. 

For subsonic flow, exit static pressure is 
specified and p, u, and v are extrapolated from 
the interior. For supersonic flow all outflow 
quantities are extrapolated. 

On the walls tangency is enforced for inviscid 
flow by extrapolating U and setting V - 0. For 
viscous flows u = v * 0. Density is extrapolated 
to the walls. Surface pressures are computed using 
the normal momentum equation: 

2 2 

(ri E + n E ) p_ + (n + n )p - - pU (n u„ n v ) 
k x s x y y 5 x y n H x \ y 5 

( 6 ) 

where U = 0 on the surface for viscous flows. 


3 = 5 x n y - 5 y n x = x y _ x y ” (cell volume) 


The contravariant velocity components U and V 
along the i- and n-grid lines are given by: 

U = E x u * 5y v ; V * n x u + n y v (5) 

S patial Differencing 

Two approaches are used to construct difference 
formulas on general curvilinear grids. In the 
finite-difference approach new Independent variables 
E, n are chosen so that 5 = 5(*,y), n = n(x,y) 
maps the original domain to a rectangle and each cur- 
vilinear coordinate to a straight line. The equa- 
tions are then transformed to ( 5 , n) coordinates 
and differenced based on even spacing A£ = An - 1. 

In the finite-volume approach the equations are 
rewritten in integral form and the divergence theo- 
rem is used to express the fluxes in terms of sur- 
face (or line) Integrals. These Integrals are then 
approximated by some integration formula. 

Using central differences these two approaches 
give rise to identical formulas for interior cells 
even for the full Navier-Stokes equations, see e.g., 
Ref. 12. Though it is not imperative, finite dif- 
ference schemes usually locate the variables at the 
nodes while finite volume schemes usually locate 
the variables at the center of the cell where they 
represent a cell-averaged quantity. The differences 
between cell-centered and node-centered schemes are 
mainly noticed at boundaries and in the transfer of 
quantities between coarse and fine grids. These 
differences will be discussed later. 


In the finite-difference code, boundary con- 
ditions are applied after all stages of the Runge- 
Kutta scheme have been completed. 

Boundary Conditions - Finite Volume Code 

For external flows the incoming Riemann 
invariant R+ = u + 2c/(y-l), the entropy, and the 
tangential velocity are specified upstream. The 
outgoing Riemann invariant R' = u - 2c/(y-l) is 
extrapolated. At a subsonic outflow the situation 
is reversed with the first three quantities 
extrapolated and R~ specified. The downstream 
total energy is evaluated assuming a constant total 
enthalpy . 

For Internal flows the total enthalpy, the 
entropy, and the tangential velocity are specified 
at the inlet and R~ is extrapolated. At the 
exit the static pressure is specified and the 
entropy, the tangential velocity, and R + are 
extrapolated . 

Variables are not defined directly on solid 
surfaces where only the normal fluxes are needed^ 
Since V = 0 on solid surfaces it follows that F 
depends only on p at the surface (Eq. 2). The 
pressure is calculated from the normal momentum 
Eq. (6). In practice a fictitious cell is placed 
inside the body. In this cell the pressure is found 
from Eq. (6), the density and normal velocity are 
set symmetric with external values, and the tan- 
gential velocity is set antisymmetric. For viscous 
flows simple extrapolation of the pressure is suf- 
ficient and both velocity components are set anti- 
symmetric. 
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Artificial Dissipation - Finite-Difference Code 


Dissipative terms consisting of fourth and 
second differences are added to prevent odd-even 
point decoupling and to allow shock capturing 
respectively. The dissipative term D in Eq. ( 1 ) 
is given by: 


Dq . m ♦ D n )q (7) 

In the finite-difference code the dissipative 
operators are non- conservative. The ^-direction 
operator is given by: 


D?q * C (V 2 q ?? - V 4 q^ 55 ) ( 8 ) 

where 



Is the arc length along the grid In the direction of 
the dissipation. The terms V 2 and V 4 are 
given by: 


v 2 = V 2 MAX J, j) 

V 4 = MAX (0, V 4 - v 2 ) 

where 




p ui..l I 2p i.j * 


p m.j + 2 p i,j + 


and 


( 10 ) 


OD 


V 2 - 0 ( 1 ) 

v 4 = O(^) (12) 


Switching functions used in Eq. (10) Increase 
V 2 slightly and switch off V 4 at shocks, 
effectively eliminating overshoots ahead of shocks. 

In smooth regions of the flow the dissipative 
terms are of third order and thus do not detract 
from the formal second-order accuracy of the 
scheme. Near shocks is large and the 

second-difference dissipation becomes locally of 
first order. 


Artificial Dissipation - Finite-Volume Code 

The artificial viscosity terms in the finite- 
volume code are similar to those in the finite- 
difference code. The dissipative term D is 
evaluated In conservation form using 


D?q = a^C (V 2 q ? - V 4 q^) (13) 

In recent versions of the code the fourth dlf- 
ference term Is computed as the second difference 
of CV 4 q^ rather than as given above. V 2 and V 4 
are given by ( 10 , 11 ) as before with typical values 
of the constants being: 

1 4- 1 

v 2 = 4 to 2 
v 4 = 128 t0 64 


Now C is given by the sum of^the spectral 
radii of the Jacoblans of E and F. Define 
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C = p( A) * P(B) = 3^ 


(15) 


An alternative is to use C = p(A) in the 
^-direction and C = p(B) in the n-direction. 

In this case the metric quantities in C represent 
a length scale normal to the direction of the dis- 
sipation, which is exactly opposite to that used In 
the finite-difference code Eq. (9). 


Multi-Stage Runge-Kutta Schemes 

An explicit Runge-Kutta scheme is used to 
advance the Euler or thin-layer Navier-Stokes equa- 
tions in time from an initial guess to a steady 
state. Given the residual R from a finite- 
difference or finite-volume representation of the 
flow equations, a k-stage scheme can be written as 

q(o) = qH 

q( 1 ) = q(o) _ a -, a t Rq(°) 


q(k) = q(o) _ a Rq(k-l) 

qn +1 . q(k) ( 16 ) 

A k-stage scheme used with central differencing 
of the inviscid equations can be made stable for a 
Courant number up to x* = k-1 , depending on the 
choice of . The values of and X* used 

here are given In Table I. For consistency we must 
have afc = 1 . For nonlinear problems the schemes 
are second-order accurate in time if <*k_i ^ 1 / 2 . 
Schemes of the form of Eq. (16) cannot be more than 
second-order accurate for any values c^. 

For efficiency both the physical and artificial 
dissipation terms are calculated once based on q(o) f 
then held constant for the remaining stages. For a 
five-stage scheme these terms are usually reevalu- 
ated after the first stage. In either case the 
multi-stage schemes have the desirable property 
that If the solution converges, i.e., Rq(o) = 0, 
then q(i) = q(o) and qn*l = qn, Independent of the 
time step. 

Spatial ly-Varlable Time Step 

Both the finite-difference and f Ini te-vol ume 
codes use a spatial ly- variable time step to 
accelerate convergence. In the finite-difference 
code the time step is given by: 
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(17) 


(1 - *5 6 «)(1 ' Vnn )R = R 


( 20 ) 
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where 


where 6^ and 6 
operators and e£ 


are standard second-difference 
and e n are smoothing parameters. 


dX z: U x | * |n X l; d y = I5yl + hyl 

This Is calculated based on the Initial con- 
ditions and is stored and not updated during the 
solution. 

In the finite- volume code the two dimensional 
time step Is given by a conservative estimate as: 

"Vj = f (1 

where 


|U| ♦ C + «y 

Aty = X (19) 

IVI * C VnJ f ny 

Multi-Stage Runge-Kutta Results 

A model problem of inviscid, transonic flow 
over a 10 percent thick circular-arc bump in a 
channel was used to investigate the effects of num- 
ber of stages and Courant number on convergence 
rate. The grid for this problem has 129 by 33 
points and is shown at the top of Fig. 1. The 
Inlet Mach number 0.675 is the one-dimensional 
choking Mach number for this case. Mach contours 
at the bottom of Fig. 1 show that although a super- 
sonic bubble and shock do develop, the two- 
dimensional flow does not choke. 

Convergence histories for the finite-difference 
scheme without multigrid in Fig. 2 show the log of 
the RMS residual in density versus iteration for two 
to five-stage schemes. Each scheme was run near its 
maximum Courant number. As expected, convergence 
rate improves with Courant number. The two- stage 
scheme required 29.6 sec of cpu time on the Cray 
X-MP 2/4 at NASA Lewis Research Center, and the 
number of iterations reached in that time is marked 
on the other curves. The higher-level schemes are 
more efficient since they obtain higher Courant 
numbers per number of stages while requiring fewer 
evaluations of the artificial dissipation. 

If the convergence criterion Is taken as a 
three decade drop in the residuals, only the five- 
stage results are fully converged at 2000 
iterations. The four-stage results, however, are 
converged to plotting accuracy. 

I mplicit Residual Smoothing 

Residual smoothing was introduced by Lerat (see 
for example Ref. 13) for use with the Lax-Wendroff 
scheme and was later applied to Runge-Kutta schemes 
by Jameson. 14 The technique involves replacing 
the residual calculated from Eq. (1) with a value 
smoothed by an implicit filter, e.g., 


Residual smoothing increases the support of the 
finite-difference scheme and thereby increases the 
stability range of the time-stepping scheme. Linear 
analysis has shown that the Runge-Kutta scheme with 
implicit smoothing may be made unconditionally sta- 
ble if e is sufficiently large. In one dimension 



gives unconditional stability where X* is the 
Courant number of the unsmoothed scheme and X is a 
larger operating Courant number. In two dimensions 
c may be reduced substantially and different values 
may be used in each direction. Turkell 5 has shown 
that since large values of e can decrease the 
convergence rate of the scheme, the best strategy 
is not to maximize the Courant number but simply to 
increase the Courant number of the unsmoothed scheme 
by a factor of two to three. 

Implicit residual smoothing requires solution 
of a scalar tridiagonal equation for each variable 
in each direction, and adds 10 to 15 percent to the 
total CPU time for a solution on a Cray X-MP. In 
the finite- volume code it is applied after each 
Runge-Kutta stage. In the finite-difference code 
it is applied after every-other stage with somewhat 
larger values of e. 

Figure 3 shows the effects of implicit smooth- 
ing on convergence rate of the five-stage scheme 
for the problem of Fig. 1. The top curve is for 
the unsmoothed scheme at a Courant number of three 
and is duplicated from Fig. 2. This scheme con- 
verges with a spectral radius of about 0.997. The 
middle curve is for a Courant number of six with 
smoothing after each stage. It converges about 
twice as fast as the unsmoothed scheme with a spec- 
tral radius of 0.994. The bottom curve for a Cour- 
ant number of nine shows a spectral radius of 
0.991. Again the convergence rate improves with 
Courant number but the amount of improvement 
decreases as the Courant number becomes large. 

Circles on the curves can be used to compare 
convergence level at a given CPU time to that of 
the two- stage scheme shown in Fig. 2. 

Multigrid Methods 

Multigrid algorithms were originally developed 
for elliptic problems. BrandtH has described 
some of the earlier versions of these algorithms. 
This approach supposes that one has an iteration 
procedure which quickly reduces the high frequency 
error but then slows down in reducing the low fre- 
quency error. One applies the basic iteration pro- 
cedure a few times to remove the high frequency 
component of the error. The remaining error is then 
fed to a coarser grid which can represent the error 
since it no longer contains high frequencies that 
would be aliased. On this coarser grid the basic 
smoothing algorithm is again used and the process 
is repeated. 
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N1 * s Method 

Nl's method Is basically a one-step Lax- 
Fredrichs scheme applied to the residual on a 
coarse grid. N1 used It to accelerate his own 
fine-grid Euler schemeJ Johnson adapted it to 
MacCormack's fine-grid scheme2 and to viscous 
flows by demonstrating that viscous terms may be 
neglected on the coarse grids. Johnson also sim- 
plified Nl's scheme considerably by replacing the 
flux Jacobian terms with temporal differences of 
the fluxes. 3 This "flux-based" coarse-grid 
scheme is outlined below. 

Given a fine-grid solution q and a fine-grid 
change 

Aq n*l = q n+l _ q n 

a coarse-grid change can be estimated using a 
Taylor series: 

iqn +2 = aqfH -1 + At ( Aq n+l) t . o(at 2 ) 

The Euler e q uat1ons are used to replace the third 
term 

(Aq"* 1 ),. = -At(E ? . F^) 


are then filled in by bilinear Interpolation and the 
fine grid Is updated using q n+2 - q n+l + Aqn+2. This 
leaves us back on the fine grid with new values of 
q and Aq, and the process can be repeated for any 
other coarse grid. In practice, the grids are 
advanced from fine to coarse with one Iteration on 
each. 

Since the right-hand side of Eq. (23) depends 
strictly on fine-grid changes. If the fine-grid 
converges the coarse-grid scheme given by (23) can- 
not change the solution. Thus, the delta form of 
Nl-type schemes retains fine-grid accuracy. For the 
same reason, all physical viscous dissipation terms 
may be neglected on the coarse-grlds without 
affecting a viscous solution on the fine grid. 

Equation (23) may be interpreted as a Lax- 
Fredrichs scheme because of the averaging of Aq 
and of the fluxes. This averaging stabilizes the 
coarse-grid scheme without the use of artificial 
viscosity. The scheme Is stable to a Courant number 
of one, which restricts the fine-grid scheme to a 
Courant number of one. Thus Nl's scheme Is attrac- 
tive for accelerating explicit HacCormack- type 
schemes but is limited to a two-stage Runge-Kutta 
scheme. In the next section we discuss a possible 
multistage Ni-type scheme. 


and the order of differentiation Is reversed 
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where 


~ ~n+l "n 
AE - E - E , etc. 


Equation (22) Is Implemented on a coarse grid with 
spacing 2hA£ and 2hAn and time step hAt, h 
* 1, 2, 4, 8 ... using: 


n+2 iff 

2h‘ML 
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] 1 +h , j 


+h 


+ Aq - At (AE - AF) 

L 1,3 J Hh.j- 
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+ Aq - At, .(-AE - AF) > 
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Finally, Eq. (22) shows that the coarse-grid 
change Aqn+ 2 W 111 p e 0 f the order of the fine- 
grid change Aq nf1 plus a small correction of 
order At. Thus each coarse-grid step gives effec- 
tively one fine grid At but with fewer opera- 
tions. This is demonstrated In Fig. 4 where con- 
vergence rates for the transonic model problem are 
compared for different numbers of grids. The fine- 
grid scheme (1 grid) is a two-stage Runge-Kutta 
scheme with a Courant number of 0.95.- Convergence 
rates for two, three, and four grids are nearly 
identical to the convergence rates for Courant num- 
bers of two, three, and four shown In Fig. 2, but 
require less CPU time. The four-grid curve In 
Fig. 4 Is about the best that we have done with 
Nl's scheme, but It Is not as fast or easy to pro- 
gram as the five-stage scheme with Implicit 
residual smoothing shown In Fig. 3. 

Multistage N1 Scheme 

The major drawback of Nl-type schemes is that 
the Lax- Fredrichs coarse-grid scheme has a Courant 
limit of one. Hence, there Is no major gain In 
using a multistage Runge-Kutta method on the fine 
grid with a large Courant number when this cannot 
be maintained on coarser grids. An alternative Is 
to use a multistage Runge-Kutta algorithm on all 
grids. This necessitates rewriting the Runge-Kutta 
scheme so that only Aq appears. Consider the 
equation 


q t = F x 


A typical stage of a Runge-Kutta scheme Is 


The fine-grid Runge-Kutta scheme Is used to 
advance the solution one time step, giving qn+l 
and Aq r,fl at each grid point. Choosing a value 
of h determines a coarser grid, and Eq. (23) is 
then used to determine Aqn+2 at those coarse-grid 
points. On the boundaries Aq n *1 is taken to be 
zero. Old values of Aq n *1 may be overwritten so 
that no additional coarse-grid storage is needed. 
Values of Aq n+2 at intermediate fine-grid points 
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(24) 

We wish to rewrite the right-hand side of (24) so 
that only Aq appears. Note that 
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Using a Taylor series expansion for F(q" + Aq j k ^ ) 
we find that an alternative to Eq. (24) Is 
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Equation (28) may now be solved on the 2h-grid 
using the Runge-Kutta scheme: 
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(29) 


and the corrections may be Interpolated back to 
the h-grid using: 


Aq 


(k+1 ) = !k±l Aq (l) 
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(30) 


- (25) 

Here only Aq’s appear, as in NT's original scheme. 
Furthermore, for a linear problem we have the stan- 
dard Runge-Kutta scheme and so the stability con- 
ditions are unchanged. 

The multistage coarse-grld scheme given by 
Eq. (25) is thus proposed as an alternative to 
Nl-type schemes that would permit larger Courant 
numbers on the coarse grids. 

Full Approximation Storage (FAS) Multlgrid Method 


On more than two grids, Eqs. (29) and (30) may 
be applied recursively. 

It can be shown that Eqs. (29) and (30) are the 
2h-grid approximation to linearizing Eq. (26) around 
the current solution qf 2 h » solving Eq. (26) for 
the error e 2h = q 2h - qf 2h , and then updating 

the current solution by q h «- q h +■ 1 2 h e 2h ‘ Th1s 

approximation Is reasonable only If the error e 2 h 
is smooth and therefore visible on the 2h-grid. 
Thus, it Is desirable to choose the parameters of 
the Runge-Kutta scheme to insure smoothing of high- 
wavenumber errors. 


Brandt’s FAS multigrid has been developed as a 
general strategy for accelerating Iterative schemes 
and may be applied directly to the multi-stage 
schemes. It results in a coarse-grid equation that 
has the same form as the fine-grid equation with 
the addition of a forcing function. The multistage 
scheme may thus be used directly on the coarser 
grids, with successively larger time steps on suc- 
cessively coarser grids. This is the primary 
advantage of Brandt’s or Jameson’s methods over 
Ni's method. Here we develop FAS multigrid using 
Brandt’s notation, then show the relationship to 
Jameson’s method. 

Brandt's notation . ^ Consider a general 
steady nonlinear equation on a grid with spacing 
parameter h. 

R h% = f h ( 26 ) 


If the fine grid converges, f h - R h q h . 0 and 
Eq. (28) becomes 

<V t * Vji, ■ vl\ • C 

or 

< q 2h> t * 0 s1nce q 2h - l l\ 

so that the coarse-grid scheme maintains fine-grid 
accuracy. Like Ni’s scheme, this is true even If 
the coarse-grid residual R 2h is different from 
the fine-grid residual R^, so that on the coarse 
grids a simple first-order artificial viscosity may 
be used and the physical viscosity may be neglected 
altogether. 


where the forcing function fh may be zero. 
Consider an approximate solution qh, evaluate 
Rh’qh and subtract it from each side of 
Eq. (26) to get the residual equation. 

R h% ~ Rh^h - fh " Rh^h 


J ameson's notation 9 . Jameson starts with a 
general unsteady equation written in f i ni te- vol ume 
form: 

d 

“ (Vq) + Rq = 0; \l = volume (31) 


Ihis may be approximated on a coarser 2h grid 
using 

R 2h q 2h * R 2h X h q h f J h (f h ‘ R hV (2 

where I^ h means interpolation from the h-grid to 
the 2h-grid. The unsteady terms may be added to 
Eq. (27) to give: 


and defines a coarse-grid forcing function as the 
difference between the coarse-grid residual and the 
sum of the fine-grid residuals: 


f 


(o) 

2h 



R 


( o) 

2h q 2h 


On the coarse grids the multistage scheme is 
implemented using 


(32) 


(q ) ¥ R q = f 

V4 2h ; t 2h 7h 2h 


(28) 


1 (1) 


= q 


(o) 

2h 


*1 at ( R 2h q 2h 


- f 


(o)\ 

2h j 


where 


or 
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At 


„<D 

q 2h 


- a (o) 

" q 2h 


R ?h<l 


2h 4 2h 


d-1) 


+ Z)h q £ 


(o) 


R a (o) 
" K 2h q 2h 


(33) 


Jameson's coarse-grid scheme Eq. (33) may be com- 
pared directly with the one based on Brandt 
Eq. (28), which Is repeated below for reference. 


n (h) 

q 2h 



wy>°\ 


(28) 


the term 1^ q^ as described above, subtract It 

from q to get the coarse-grid correction, then 
2h 

transfer the correction to the fine grid using 
bilinear interpolation. If more than one 
Runge-Kutta sweep is performed on a coarse grid then 

I^ h q h must be stored or reconstructed. 

Thus the main difference between the finite- 
difference and f ini te-volume multigrid schemes is 
the restriction operators which have been chosen 
appropriately for each scheme. We suspect that the 
differences have minor effects on convergence 
rates, but have not quantified those effects. 


It is evident that the two schemes are equiva- 
lent in form but differ in the interpolation or 
grid-transfer operators. 

Grid-Transfer Operators 

Before discussing these operators we note some 
differences in coarse-grid data locations between 
node-centered (finite-difference or F.D.) and cell- 
centered (finite-volume or F.V.) schemes. In 
node-centered schemes, coarse-grid nodes coincide 
with fine-grid nodes (Fig. 5(a)) so that variables 
can be transferred between grids directly without 
interpolation. This process is termed "injection." 
In cell-centered schemes, however, coarse-grid 
cells are made up of several fine-grid cells, and 
the centers of the coarseand fine-grid cells are 
not coincident (Fig. 5(b)). Thus, it is always 
necessary to use some form of interpolation between 
grids with cell-centered schemes. The grid transfer 
operators used in the FD and FV codes are listed 
below: 

Restriction operators . For the term I^ h q h the 

the FD code uses injection and the FV code uses a 
volume-weighted average. 

FD: P 2 h * Ph 


Coarse-Grid Boundary Conditions 

Most of the boundary conditions described ear- 
lier Involve extrapolations or one-sided dif- 
ferences and do not maintain fine-grid accuracy 
when applied directly on coarse grids. It is pos- 
sible to define a boundary condition forcing func- 
tion similar to Jameson's Interior forcing function 
Eq. (32), but this increases programming com- 
plexity. Instead, coarse-grid boundary conditions 
are computed with coarse-grid accuracy, but only 
the change in boundary values during the Runge- 
Kutta cycle is transferred back to the fine grid, 
i .e. , 

_ h . (k) (o) . 

q hBC * q hBC + *2h q 2hBC ' q 2hBC ) (36) 

This formulation maintains consistency in that if 
the interior scheme converges the coarse-grid scheme 
cannot change the fine-grid boundary values. 

In the finite-difference code, boundary con- 
ditions are updated after every complete Runge- 
Kutta cycle. In the original version of FL052MG, 
boundary conditions were frozen on coarse grids. 

In the present version of FL052MG, boundary con- 
ditions are updated after every stage of the Runge- 
Kutta scheme and after every grid transfer. 


FV: 


q 2h 


yVh 

V 2h 


(34) 


2h For the restriction of the fine-grid residual, 
I^ n (R h q h - f^)» Brandt recommends a weighted-average 

over several grid points. Although the FD code 
could use injection here we have obtained much 
better convergence rates using an unweighted average 
over nine neighboring nodes. The FV code uses a 
sum over the four fine-grid cells that make up the 


coarse-grid 

2h 

cell. 

4 



FO: I, 

< R h q h - V 

= 9 

L 

9 nodes 

< R h q h 

"= >? 

< R h q h - V 

= 

Z 

4 cells 

< R h q h 


Programming Considerat i ons 

FAS multigrid is usually programmed by storing 
all the grids end-to-end in one long, singly- 
dimensioned array. This Increases the required 
storage by about 4/3 in two-dimension or 8/7 In 
three-dimension. Brandt typically accesses these 
long arrays by indirect addressing which can be 
difficult to program and vectorize. Jameson over- 
comes these difficulties in FL052MG by defining 
singly-dimensioned arrays In the main program but 
working with multiply-dimensioned arrays In the 
subroutines which have the array size and starting 
location passed through their argument list. This 
technique allows existing nonmultigrid subroutines 
to be converted to multigrid with minimal changes. 
Nevertheless, FAS multigrid is considerably more 
difficult to program than Ni's method. 

Full Multigrid 


(35) 

P rolongation operators . For the expression 
q h *■ q h * 1 2h (q 2h “ 1 h hq h ) * b0th COdeS calculate 


Full multigrid ( FMG) combines successive grid 
refinement with FAS multigrid. The solution is 
started on the coarsest grid and is Iterated a few 
times using the fine-grid scheme, possibly with FAS 
multigrid as well. The solution is next Inter- 
polated to the next-finer grid where it provides a 
good initial guess, then the process is repeated 
until the finest grid is reached and the 
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solution converges. Experience has shown that It 
Is wasteful to drive the residual to zero on coarse 
FMG grids. We have found that 25 Runge-Kutta cycles 
are sufficient, with more Iterations used only on 
the finest grid. 

Multlgrld Cycle 

In the finite-difference code we use a V or 
sawtooth- eye le with one Runge-Kutta Iteration on 
each grid. A typical cycle Is diagrammed In Fig. 6. 
Turkel's experience has shown some Increase In 
efficiency by doing one Interation on the h-grld, 
two iterations on the 2h-grid, and three Iterations 
on all coarser grids. More Iterations on the 
coarsest grid may actually decrease the convergence 
rate. 

Multlgrld Results 

Figure 7 shows FMG convergence rates for the 
problem of Fig. 1 using the finite-difference 
code. Here 300 iterations were run on each of four 
successive grids to show the convergence rate on 
each. Normally only 25 iterations would be run on 
the first three grids. The coarsest grid (17 by 5) 
was run with the four-stage scheme alone and has a 
spectral radius of 0.962. The second and third 
grids have spectral radii of 0.937, and the finest 
grid (129 by 33) converges with a spectral radius 
of 0.95. 

Figure 8 shows the FMG convergence rate of 
Turkel's version of Jameson's f ini te- vol ume code 
FL052MG for external flow about an RAE2822 super- 
critical airfoil with M^ - 0.75 and a = 3°. 

A C-type grid was used with 160 by 32 points on the 
finest grid. Thirty Iterations were run on both 
the coarse and medium grids, with two and three 
grid levels used respectively. Ninety Iterations 
were run on the finest grid with four grid levels, 
giving a spectral radius of about 0.763. This is a 
typical convergence rate for this code for inviscid 
external flows on uniform grids. 

Figure 9 compares convergence rates for Ni 1 s 
method, the finite-difference Runge-Kutta scheme 
alone, and the finite-difference and finite- volume 
Runge-Kutta schemes with FAS multigrid. The Ni 
multigrid results were run with four grids at a 
Courant number of 0.95, and show a spectral radius 
of 0.995. The four-stage Runge-Kutta results were 
run with Implicit residual smoothing at a Courant 
number of 5.2 and show a spectral radius of 0.994. 
The finite-difference FAS multigrid results were 
run with the same Runge-Kutta parameters as above 
but with 25 FMG iterations on each coarse grid and 
300 FAS iterations on the finest grid. The spectral 
radius of the FAS multigrid scheme is 0.943. The 
finite- volume FAS multigrid results were run with a 
four-stage scheme at a Courant number of 6.0 and 
also with a five-stage scheme at a Courant number of 
7.5. These results show spectral radii of 0.855 and 
0.822 respectively. The differences in convergence 
rate between the finite-difference and finite- 
volume codes are probably due to the additional 
enthalpy damping step used in the finite- volume 
code . 

Circles on Fig. 9 show equal CPU times of 
16.7 sec on the Cray X-MP. In terms of CPU time 
required to reach a certain convergence level, the 
implicitly-smoothed Runge-Kutta scheme alone is 
marginally faster than Ni's scheme. The two FAS 


multigrid schemes are similar to each other in per- 
formance and are approximately 8.5 times faster 
than Ni's scheme or the implicitly-smoothed scheme 
alone. 

Error Smoothing Versus Time-Marching 

In the original work on multlgrld for elliptic 
equations, point Jacobi or Gauss-Seidel iteration 
schemes were used because of their ability to smooth 
high-frequency errors. The fact that these schemes 
could be interpreted as time-marching schemes was 
considered irrelevant. In contrast Ni's scheme has 
been referred to as "hyperbolic multigrid" with its 
main purpose being to advance rapidly in time using 
coarse-grid information. 

It is not clear which interpretation is more 
appropriate for analysis of FAS multigrid schemes 
applied to hyperbolic problems. Jameson 9 bases 
his analysis on the smoothing properties of the 
Runge-Kutta scheme. On the other hand, 

Jesperson 16 has shown that multigrid solutions 
possess time-accurate properties. Thus FAS may 
work because it allows larger time-steps on the 
coarser grids. 

Using a stability analysis of the linear one- 
dimensional convection equation with fourth- 
difference artificial dissipation it is possible to 
choose the multistage scheme parameters aj<, y, 
e, and \ to maximize the smoothing and the time 
step. Jameson has published several sets of these 
parameters in Ref. 9. Ouring the current work 
Schaffer developed an optimization code that chooses 
these parameters to minimize the area under the 
ampl i f ication factor curve. 

Experience with these "optimal parameter 
schemes" has been inconclusive. It is clear that 
both large Courant numbers and high smoothing lead 
to fast convergence. However, schemes with equal 
Courant numbers but very different one-dimensional 
smoothing properties often have similar convergence 
rates. We suspect that the effects of stretched 
grids on the artificial viscosity and the effects 
of applying implicit residual smoothing as a 
sequence of one- dimensional operators are such that 
the ampl if ication factor for the two dimension codes 
do not look much like the one- dimensional model 
results. Also, Jameson has shown that amplification 
factors can vary considerably from grid to grid!7 
making it difficult to predict the behavior of the 
overall multigrid scheme. 

Viscous Results 

Experimental data for the interaction of an 
oblique shock wave with a laminar boundary layer 
have been published by Hakkinen et al. 18 This 
case has been computed by many researchers, notably 
by MacCormack and Baldwin 19 using a 32 by 32 mesh. 

Here we have computed this flow with a 113 by 
41 mesh. The free-stream Mach number is 2. The 
upper boundary was treated as an inviscid wall 
which was bent 3.091° to generate an oblique shock 
that intersects the lower wall at 0.16 ft at an 
angle of 32.585°. The grid has a constant x-spacing 
of 0.003 ft and a y-spacing at the wall of 
0.0001 ft, stretching geometrically through the 
boundary layer to a constant spacing above. At the 
wall the grid aspect ratio is 30:1. 
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Figure 10 shows an Euler solution computed on 
the viscous grid. Static pressure contours are 
shown at the bottom and the lower-wall pressure 
distribution Is compared to experimental data at 
the top. The Invlscld lower wall has no effect on 
the flow except to reflect the Incident shock. The 
computed pressure distribution shows some overshoot 
after the reflected shock but a good comparison 
with the overall shock strength. 

Convergence histories for the Ni, Runge-Kutta, 
and FAS multigrid schemes are shown In Fig. 11. All 
three schemes show fast Initial convergence for this 
purely supersonic flow, but all slow down abruptly 
after a three-to-four decade drop In the residuals, 
presumably because of the very fine grid at the 
wall. It appears that the multigrid schemes par- 
allel the convergence behavior of the fine-grid 
scheme but at a faster rate. The N1 scheme shows a 
surprising dip and jump In convergence near the end. 
Circles showing convergence levels at equal CPU 
times of 22.3 sec Indicate that the three schemes 
are comparable in speed for this case. 

Figure 12 shows the thin-layer Navler-Stokes 
solution for this problem. Flow conditions were 
chosen to give a Reynolds number of 2.96x105 at 
the shock (x - 0.16 ft). Static pressure contours 
at the bottom of Fig. 12 show not only the Incident 
and reflected shocks but a leading edge shock that 
reflects from the upper wall and compression waves 
generated by the separation bubble. The lower-wall 
static pressure distribution Is compared to the 
experimental data at the center of Fig. 12 and shows 
good agreement after the leading-edge and reflected 
shocks, but slightly under predicts the plateau 
pressure In the separated region. Computed skin 
friction shown at the top of Fig. 12 shows very good 
agreement with the data. Including separation and 
reattachment points. Note that skin friction was 
not measured in the reverse flow region. 

Viscous flow convergence rates shown In Fig. 13 
are slower than for the Euler case and show the same 
tendency to slow down, but not so abruptly. Circles 
at equal CPU times of 28 sec show that the Runge- 
Kutta scheme with Implicit smoothing is slightly 
faster than Ni's scheme, but that the FAS multigrid 
scheme Is about 2.1 times faster. 

Cone 1 usl on 

In this work we have developed and compared the 
Ni, Jameson, and Brandt types of multigrid- schemes. 
Each scheme uses the Runge-Kutta method as the basic 
algorithm and that method has been studied without 
multigrid as well. 

We have demonstrated that the efficiency of 
the Runge-Kutta scheme Increases with number of 
stages and that convergence rate increases with 
Courant number. We have shown that Ni's scheme 
gives one fine-grid time step on each coarse grid 
but that its Courant limit of one limits the effec- 
tiveness of the scheme. We propose a multistage 
coarse-grid scheme to improve the Courant limit. 
Using implicit residual smoothing to increase the 
Courant limit of the Runge-Kutta scheme without 
multigrid gives a scheme slightly more efficient 
than Ni's. 


with differences between Interpolation methods for 
the cell-centered and node-centered schemes. For 
an invlscld transonic model problem the finite- 
difference multigrid scheme is about 8.5 times 
faster than the best NI or Runge-Kutta scheme with- 
out multigrid. Jameson's finite-volume code Is 
about three times faster than that in terms of 
Iterations, probably due to the additional enthalpy 
damping step. The two multigrid codes are about 
equal In terms of CPU time. 

Convergence rates of all methods decrease on 
highly-stretched grids. For a shock-boundary layer 
Interaction problem the finite-difference multigrid 
code was still 2.1 times faster than the best NI or 
Runge-Kutta scheme without multigrid. We believe 
that significant Improvements can still be made In 
multigrid convergence rates for viscous flows. 
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Figure 1.- Top: 129x33 grid over a 10% thick circular arc bump. 
Bottom: Mach contours, M jn = 0.675, sonic line is dashed. 


Figure 2.- Comparison of convergence rates for k-stage schemes at 
various Courant numbers A. 
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Figure 3.- Comparison of convergence rates at Courant No. A. Figure H.- Comparison of convergence rates for N-grid Ni -Schemes. 

5-stage scheme, 2 evaluations of dissipation. Implicit 
smoothing with parameter e after each stage. 





Figure 5.- Relationship between fine- and coarse-grid 

PO j mjj 

(A) NODfc-CbNIERED (FINITE DIFFERENCE) SCHEMES. 

(b) Cell-centered (finite volume) schemes. 



25 CYC 25 CYC 25 CYC 300 CYC 

1 GRID 2 GRIDS 3 GRIDS <4 GRIDS 



Figure 7,- Fiji i multigrid convergence rates (MAX AND RMS density 
residual) four stages, A= 5.2. 


Figure 6.- Typical FMG cycle for the finite-difference 
code. 
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Figure 8.- Convergence history for inviscid flow about a RAE 2822 Figure 9 ' CoMPAR ' SON 0F convergence rates for five different schemes. 

AIRFOIL USING THE F I N I TE -VOLUME MULT I GRID METHOD. 




Figure 10.- Oblique shock reflection, Euler solution. Figure 11.- Comparison of convergence rates for inviscid flow. 

Top: Comparison of computed and measured static pressures. 

Bottom: Pressure contours. 
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